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Battery  cell  life  depends  critically  on  how  the  cell  is  used.  Therefore,  battery  chargers  and  battery  man¬ 
agement  systems  must  be  designed  to  control  cell  usage  carefully.  In  order  to  design  optimal  battery 
controls  that  effect  a  tradeoff  between  cell  performance  (in  some  sense)  and  cell  life,  a  model  of  cell 
degradation  is  necessary.  This  model  must  be  simple  and  incremental  in  order  to  be  implemented  by 
an  inexpensive  microcontroller.  This  paper  takes  a  first  step  toward  developing  such  a  controls-oriented 
comprehensive  cell  degradation  model  by  deriving  a  reduced-order  model  of  a  single  mechanism:  the 
growth  process  of  the  solid-electrolyte  interphase  (SEI)  layer,  along  with  the  resulting  resistance  rise  and 
capacity  loss.  This  reduced-order  model  approximates  a  physics-based  PDE  model  from  the  literature,  is 
simple  and  accurate,  and  may  be  used  in  optimal  strategies  for  controlling  lithium-ion  batteries. 

©  2012  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Many  applications  of  battery  packs  require  high  power  and  long 
life.  These  two  objectives  are  at  odds  with  each  other  since,  in 
general,  high  charge  and  discharge  rates  accelerate  aging.  It  is  our 
ultimate  goal  to  design  methods  to  control  a  battery  pack  to  effect 
an  optimum  tradeoff  between  performance  and  degradation,  and 
to  do  so  we  will  need  mathematical  models  of  the  dominant  cell 
degradation  mechanisms.  These  models  must  be  simple  enough 
to  be  quickly  executed  on  an  inexpensive  embedded  systems 
processor. 

Presently,  many  battery  management  systems  (BMS)  track 
macroscopic  indicators  of  aging  such  as  capacity  fade  and  power 
fade  but  few,  if  any,  track  physical  indicators  of  aging  such  as 
degree  of  solid-electrolyte-interphase  (SEI)  layer  growth  on  anode 
particles,  degree  of  lithium  plating  on  anode  particles,  degree  of 
cathode  dissolution,  and  so  forth.  However,  identical  combinations 
of  capacity  loss  and  resistance  rise  might  be  achieved  by  travel¬ 
ing  different  paths— by  exciting  different  physical  mechanisms  of 
degradation— yielding  potentially  different  resultant  optimal  con¬ 
trol  strategies  for  operating  the  battery  pack.  Therefore,  knowing 
cell  resistance  and  capacity  alone  is  not  sufficient  to  extend  cell  life 
to  the  maximum  extent  via  control. 
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We  believe  that  by  modeling  the  dominant  physical  degradation 
mechanisms  that  occur  in  a  cell,  and  then  by  using  those  mod¬ 
els  in  an  optimized  predictive  control  algorithm,  battery  life  can 
be  extended.  For  a  practical  implementation,  this  requires  that  the 
degradation  model  be  simple  and  incremental.  We  consider  a  model 
to  be  simple  if  it  is  described  by  a  small  number  of  ordinary  differ¬ 
ence  equations,  and  incremental  (or  recursive)  if  it  predicts  a  future 
degradation  state  based  on  a  present  degradation  state  and  a  pro¬ 
posed  cell  input  current.  The  model  will  necessarily  be  a  function  of 
cell  temperature,  state-of-charge  (SOC),  and  possibly  other  factors 
as  well. 

Degradation  leading  to  capacity  loss  and  resistance  rise  can 
occur  either  due  to  mechanical  stress  factors  or  chemical  side 
reactions  [1-8].  Various  schemes  have  been  introduced  to  model 
side  reactions,  including  that  of  Darling  and  Newman,  who  first 
introduced  the  concept  of  modeling  parasitic  effects  in  lithium-ion 
batteries  by  modeling  solvent  oxidation  side  reactions  [3].  More 
recently,  Ramadass  and  colleagues  have  proposed  a  model  that 
describes  the  formation  and  growth  of  an  SEI  layer  on  anode  solid 
particles  during  charging,  that  uses  solvent  reduction  as  the  main 
side  reaction  mechanism  for  degradation  [9].  The  present  paper 
builds  on  and  extends  this  work  to  propose  a  simple  incremental 
model  of  SEI  growth  and  associated  capacity  loss  and  resistance 
rise.  The  proposed  model  may  then  be  used  in  an  optimal  control 
scheme  that  will  reduce  side  reactions  and  limit  cell  degradation. 
(The  control  task  is  the  subject  of  present  research,  and  will  be 
reported  in  a  later  paper.  The  methodology  presented  in  [  1 0]  is  one 
possible  approach.)  While  various  efforts  have  been  made  to  cre¬ 
ate  reduced-order  models  of  ideal-cell  dynamics  (i.e.,  models  that 
do  not  consider  degradation),  including  single  particle  models  by 
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Nomenclature 

a  specific  surface  area  of  porous  electrode,  nrr1 

brug  Bruggeman  exponent,  unitless 
c  concentration  of  Li  or  Li+  ions,  mol  nrr3 

D  diffusion  coefficient,  m2  s-1 

F  Faraday’s  constant,  96  487  C  mol-1 

10  exchange-current  density  for  intercalation  reaction, 
Am-2 

z0  s  exchange-current  density  for  side  reaction,  A  nrr 2 
Jj  local  volumetric  current  density  for  intercalation 

reaction,  Am-3 

Js  local  volumetric  current  density  for  side  reaction, 

Am-3 

k  rate  constant  of  electrochemical  reaction, 

Am5/2  mol-3/2 
L  length  of  cell,  m 

Q.  cell  capacity,  C 

R  particle  radius,  p,m 

flfiim  film  resistance  at  the  electrode/electrolyte  interface, 
Q  m2 

Rg  universal  gas  constant,  8.314 J  mol-1  K-1 

T  temperature,  K 

t  time,  s 

£+  transport  number,  unitless 

U  local  equilibrium  potential,  V 

x  length  dimension,  m 

Greek 

aa,  oic  anodic  and  cathodic  transfer  coefficients  of  electro¬ 
chemical  reaction,  unitless 
e  volume  fraction  of  a  phase,  unitless 

0  local  potential  of  a  phase,  V 

1 1  local  overpotential  driving  electrochemical  reac¬ 
tion,  V 

k  conductivity  of  electrolyte,  S  m_  1 

0  state-of-charge  or  stoichiometry  of  electrode 

a  conductivity  of  electrode,  S  nr1 

p  density  of  active  material,  kg  nrr3 

8  film  thickness,  m 

Subscrip  t/superscrip  t 

e  pertaining  to  the  electrolyte  phase 

s  pertaining  to  the  solid  phase 

n  pertaining  to  the  negative  electrode 

p  pertaining  to  the  positive  electrode 


Chaturvedi  et  al.  [11]  and  one-dimensional  models  by  Smith  et  al. 
[  1 2  ],  we  believe  this  to  be  the  first  attempt  to  create  a  reduced-order 
model  to  describe  a  side  reaction. 

We  note  that  the  field  of  creating  mathematical  models  of  cell 
degradation  mechanisms  is  still  developing.  The  accuracy  of  the 
reduced-order  model  will  be  no  better  than  the  accuracy  of  the 
underlying  partial-differential-equation  (PDE)  model  from  which 
it  was  created.  So,  it  is  important  for  the  reader  to  understand  that 
the  purpose  of  this  paper  is  to  propose  a  methodology  for  reducing 
the  order  of  a  PDE  model  of  a  degradation  mechanism  to  be  able  to 
approximate  the  PDE  model  with  one  that  can  be  implemented  on 
a  simple  microprocessor,  such  as  in  one  found  in  a  battery  manage¬ 
ment  system.  For  purpose  of  illustration,  we  apply  the  methodology 
to  the  specific  cell  degradation  mechanism  of  SEI  layer  growth  as 
modeled  by  Ramadass  et  al.  [9].  The  methodology  could  also  be 
applied  to  other  similar  models;  for  example,  as  newer  and  more 
accurate  models  are  reported  in  the  literature. 


In  particular,  our  order-reduction  method  uses  volume  averag¬ 
ing  to  create  an  algebraic  “zero-dimensional”  or  “0-D”  model  of 
the  infinite-order  PDE  model  proposed  in  [9].  This  reduced-order 
model  (ROM)  of  the  SEI  growth  mechanism  is  a  first  step  toward 
creating  a  complete  coupled  reduced-order  model  of  all  dominant 
cell  degradation  mechanisms,  which  then  could  be  used  in  an  opti¬ 
mal  control  scheme. 

This  paper  will  proceed  by  reviewing  the  model  of  SEI  layer 
growth  proposed  in  [9].  A  reduced-order  volume-averaged  approx¬ 
imate  model  is  then  derived.  Results  of  pulsed-current  simulations 
of  the  full-order  model  and  the  0-D  reduced-order  model  are  pre¬ 
sented  and  compared,  results  are  discussed,  and  conclusions  are 
made. 

2.  Original  model 

Changes  at  the  electrode/electrolyte  interface  due  to  side  reac¬ 
tions  between  the  anode  and  electrolyte  are  considered  to  be  one 
of  the  primary  causes  of  anode  aging  [8].  There  are  a  large  num¬ 
ber  of  reduction  reactions  that  can  lead  to  the  deposition  of  solid 
products  on  the  anode,  and  these  are  less  well  understood,  being 
very  dependent  upon  the  composition  of  the  electrolyte  solution 
[13].  Thus  the  assumption  for  the  original  Ramadass  model  is  that 
the  side  reaction  is  considered  to  be  a  consumption  of  the  sol¬ 
vent  species  and  lithium  ions,  which  will  form  compounds  such  as 
Li2C03,  LiF,  Li20,  and  so  forth,  depending  on  the  nature  of  the  sol¬ 
vent.  Previous  studies  of  the  SEI  layer  formation  on  lithiated  carbon 
have  shown  that  there  is  probably  significant  porosity  in  the  film, 
and  this  makes  it  reasonable  to  assume  that  the  SEI  layer  contin¬ 
ues  to  grow  as  the  solvent  diffuses  through  the  layer  during  charge 
[  1 4-1 7  ].  The  assumption  of  the  ongoing  formation  of  the  SEI  layer  is 
also  supported  by  the  research  of  Aurbach  and  colleagues  [14],  who 
propose  that  the  intercalation  of  lithium  into  the  graphite  anode 
leads  to  increase  in  the  lattice  volume,  which  in  turn  stretches  the 
SEI  layer,  causing  it  to  fracture  and  to  expose  more  of  the  anode 
to  the  electrolyte,  fueling  the  side  reaction,  and  contributing  to  SEI 
formation. 

As  our  goal  is  to  create  a  high-fidelity  reduced-order  model  of 
the  degradation  model  in  [9],  we  adopt  the  same  assumptions  as 
they,  which  were: 

1.  The  main  side  reaction  is  due  to  the  reduction  of  an  organic 
solvent,  expressed  as  S  +  2Li+  +  2e_  ->►  P,  where  “S”  refers  to  the 
solvent  and  “P”  to  the  product  formed  in  the  side  reaction. 

2.  The  reaction  occurs  only  during  charging  of  the  cell. 

3.  The  products  formed  are  a  mixture  of  different  species,  result¬ 
ing  in  averaged  mass  and  density  constants  used  in  the  later 
equation  describing  the  formation  and  growth  of  the  SEI  film. 

4.  The  side  reaction  is  assumed  to  be  irreversible  and  Urefs  is  chosen 
to  be  0.4  V  versus  Li/Li+. 

5.  The  initial  resistance  of  the  SEI  layer  developed  during  cell  for¬ 
mation  is  100  £2  cm2. 

6.  There  is  no  overcharge  reaction  considered  (i.e.,  lithium  plating 
is  not  modeled). 

We  have  somewhat  relaxed  Assumption  2,  allowing  side  reac¬ 
tions  to  occur  during  rest  intervals  also  (and  even  during  discharge, 
which  we  do  not  report  here).  The  reduced-order  modeling  method 
that  we  present  easily  incorporates  either  variant  of  the  assump¬ 
tion. 

The  SEI  growth  model  of  [9]  is  tightly  coupled  with  a 
Doyle-Fuller-Newman  porous-electrode  style  model  of  ideal-cell 
dynamics  [4],  which  assumes  that  the  solid  and  electrolyte  phases 
are  considered  continuous  and  gives  no  consideration  to  the  under¬ 
lying  microstructure  of  the  electrode.  As  porous-electrode  theory 
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is  well  known,  but  as  there  are  variants  in  its  implementation,  we 
present  those  equations  in  Appendix  A,  and  focus  only  on  the  SEI 
growth  model  of  [9]  here. 

For  the  negative  electrode  the  local  volumetric  transfer  current 
density  Jtotal  is  given  by  a  sum  of  the  intercalation  current  density 
Ji  and  the  side  reaction  current  density  Js , 


./total  —Jl+Js, 


(1) 


where  Ji  is  computed  via  the  Butler-Volmer  electrochemical  kinetic 
expression 

nF  \ 

—  rinj  -  exp 
which  is  driven  by  the  overpotential 

Vn  =  4>s  -<Pe-  Uf  -  (3) 

Un 

where  i0,n  is  the  exchange  current  density  and  U^ef  is  the  equilib¬ 
rium  potential  which  is  evaluated  as  a  function  of  the  solid  phase 
concentration  at  the  surface  of  the  particle. 

The  kinetics  of  the  side  reaction  are  described  using  a  Tafel  equa¬ 
tion,  which  assume  that  the  side  reaction  is  considered  irreversible, 


Otc,nF 

'  V7' 


On 


(2) 


Ji  =  <3ni0,n 


exp 


Oin, 


Js  =  -io,sa„  exp  (  -^frh 


and  the  side  reaction  overpotential  is  described  as 

ris  =  <Ps-</>e-Urels-J-^Rmm. 

un 


(4) 


(5) 


alone,  neglecting  local  variations  in  electrolyte  and  solid  sur¬ 
face  concentration.  The  estimated  value  of  Js  then  corresponds 
to  a  suddenly  applied  current  pulse  of  magnitude  iapp,  which  is 
constant  over  some  time  interval  At. 

2.  The  intercalation  current  density  and  the  side-reaction  current 
densities  are  uniform  over  the  anode.  This  allows  us  to  state  that 
the  total  reaction  current  density  Jtotai  is  related  to  the  applied 
cell  current  iapp  by  the  following  relationship: 

Jtotal  =  -^,  O) 

where  the  volume  of  the  active  material  is  described  by 
Vol„  =L„A. 

3.  The  anodic  and  cathodic  charge-transfer  coefficients  are  equal 
(a  a  =ac  =  0.5). 


From  the  above  assumptions,  an  incremental  degradation  model 
can  be  formulated  as  follows.  First,  at  any  point  in  time,  the  lithia- 
tion  state  of  the  anode  is  calculated  as 

@n  —  @n,  min  +  SOCcell  (0„  ,max  @n,  min  ),  (10) 

where  0n, max  and  0nmin  are  the  stoichiometric  limits  of  anode  lithi- 
ation  (i.e.,  the  value  of  0  in  Li0C6  when  the  cell  is  fully  charged 
and  discharged,  respectively).  SOCcell  is  a  value  between  zero  and 
one,  which  indicates  the  cell  state-of-charge.  Then  U„ef  is  calculated 
from  0n  using,  for  example,  Eq.  (15). 

We  will  ultimately  iterate  to  find  Js.  We  can  initialize  its  value 
to  zero  and  calculate  the  intercalation  current  density  from  Eqs.  (1 ) 
and  (9). 


Once  the  side  reaction  current,  Js,  has  been  calculated,  the 
change  in  the  film  thickness  during  charging  can  be  calculated 
by 


d^fiim  =  Mp 
3 1  anpPF 


(6) 


where  MP  is  the  average  molecular  weight  of  the  constituent 
compounds  of  the  SEI  layer  and  pP  is  the  average  density  of  the 
constituent  compounds.  This  allows  the  overall  film  resistance  to 
be  calculated  as 

Rfilm  =  RsE.  +  ^,  (7) 

Kp 


where  Rs ei  is  the  initial  film  resistance  that  is  produced  during  the 
formation  period  of  the  battery,  and  kp  is  the  conductivity  of  the 
film. 

In  addition  to  the  resistance  change,  there  is  a  capacity  loss 
caused  by  the  side  reaction  current  during  charge,  leading  to  capac¬ 
ity  changing  via  the  relationship 

Cln 

JsAdx.  (8) 


9Q 

3 1 


-s 


Ji 


-*app 

Vol„ 


-Js. 


From 7/  and  assumption  3,  we  can  manipulate  Eq.  (2)  to  find  the 
overpotential  for  the  intercalation  reaction  to  be 


On  = 


2  RgT 
F 


asinh 


Then,  from  Eqs.  (3)  and  (5),  the  overpotential  for  the  side  reaction 
is: 


Os  —  <fis  -  <fie  -  UKf  s  - 

=  rin  +  Urnef  -UTeftS. 


./total  ^fil  m 
On 


We  note  that  the  film  resistance  cancels  from  the  calculation.  We 
can  now  calculate  an  updated  estimate  of  the  side-reaction  current 
density  as 


Js 


-io,san  exp 


In  total,  we  have  the  reduced-order  model 


3.  Simplifying  the  model 

To  effect  an  optimal  control  strategy,  the  processors  in  battery 
chargers  and  battery  management  systems  must  be  able  to  calcu¬ 
late  the  side  reaction  current  density Js  very  quickly  and  accurately. 
Solving  the  coupled  PDE  equations  described  above  is  too  compli¬ 
cated  for  such  a  process.  The  Js  model  needs  to  be  much  faster  and 
simpler.  In  this  section,  we  present  a  simpler  incremental  model 
for  calculating  Js,  Kfilm,  and  Q. 

To  create  a  volume-averaged  0-D  reduced-order  model,  three 
additional  assumptions  are  made: 

1.  The  cell  is  always  in  a  quasi-equilibrium  state,  allowing  the 
exchange  current  density  i0,n  to  be  calculated  from  the  cell  SOC 


Js  =  -io,san  exp 

+uf-urel 


-F 


2RgT  •  u 
— -f—  asinh 


=  -lo, sfln  exp 


— *'app/Vol  n  —  Js 

2anio 

-F  (U/,ef  —  Uref  s 
2  RgT 


Xe-P(-aSinh(-',2/0V*-'1'))’ 

which  we  iterate  until  convergence  (we  find  that  fewer  than  ten 
iterations  are  generally  necessary— the  results  presented  in  this 
paper  use  exactly  ten  iterations). 

Once  we  have  solved  for  Js  it  can  then  be  incorporated  into  incre¬ 
mental  equations  for  film  resistance  and  capacity  loss.  It  is  assumed 
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that Js  is  constant  over  some  small  time  interval  At,  and  is  denoted 
as  JS[N]  for  the  Nth  interval.  Since,  according  to  Eq.  (6),  the  change  in 
film  thickness  is  proportional  to  Js ,  we  can  arrive  at  an  incremental 
equation  for  film  thickness  as 

SrnmW  =  Sam[N  -  1  ]  -  -  1  ], 

anppt 

where  At  is  the  equation  update  period,  and  noting  that  the  sign  of 
Js  is  negative.  This  result  can  be  used  to  calculate  the  film  resistance 
as 

iWJV]  =  iWN  -  1  ]  -  -  1  ]• 

Similarly,  we  can  discretize  the  capacity  equation  (8)  to  find  that 
Q[N]  =  Q[N-l]  +  (iU„At)/s[N-l]. 


In  summary,  the  proposed  reduced-order  model  (ROM)  equa¬ 
tions  are: 

@n  —  $n, min  +  SOCcep  max  —  $n, min)  (1 1 ) 


Table  1 

Electrode  parameters  for  simulation. 


Symbol 

Units 

Anode 

Separator 

Cathode 

l 

[Jim 

88 

20 

80 

R 

[Jim 

2 

- 

2 

A 

m2 

0.0596 

0.0596 

0.0596 

a 

Sm-1 

100 

- 

100 

Ss 

- 

0.49 

- 

0.59 

se 

- 

0.485 

1 

0.385 

brug 

- 

4 

- 

4 

csmax 

mol  m3 

30  555 

- 

51555 

c° 

mol  m3 

1000 

1000 

1000 

$i,min 

- 

0.03 

- 

0.95 

$i,max 

- 

0.886 

- 

0.487 

Ds 

m2  s_1 

3.9  x  10“14 

- 

1.0  x  10“14 

De 

m2  s_1 

7.5  xlO-10 

7.5  xlO”10 

7.5  xlO”10 

t? 

- 

0.363 

0.363 

0.363 

k 

Am5/2  mol_3/2 

4.854  x  10“6 

- 

2.252  x  10“6 

OLa 

- 

0.5 

- 

0.5 

OLc 

- 

0.5 

- 

0.5 

t^ref ,s 

V 

0.4 

- 

- 

io,s 

Am-2 

1.5  x  10“6 

- 

- 

Js[N]  =  -i0?san  exp 


'-F  (^ef-Uref,s 

2  RgT 


x  exp  (  -asinh 


-iaPP[N]/Voln  -UN] 
2ani'o 


(12) 


RfilmlN]  =  Rfil m[N  -  1]  -  J^~Js[N  -  1] 
anpptKp 

QIN]  =  QIN  -  1  ]  +  {ALnAt)Js[N  -  1  ]. 


(13) 

(14) 


4.  Comparing  the  models 

The  validity  of  this  reduced-order  model  depends  first  on  the 
accuracy  of  the  underlying  partial  differential  equation  model, 
which  we  assume  here  to  be  exact.  It  then  depends  on  how 
closely  the  reduced-order  approximation  of  Js  matches  the  exact 
calculation  of  Js.  In  this  section,  results  from  both  the  full  and 
reduced-order  models  for  Js  are  compared.  In  addition,  the  error 
in  Js  is  tracked  over  time  in  order  to  see  whether  the  assumption 
that  Js  is  constant  over  a  small  period  of  time  is  reasonable. 

To  compare  the  PDE  and  reduced-order  models,  we  conducted  a 
series  of  simulations.  In  each  simulation,  the  cell  was  initially  at  rest. 
A  sudden  pulse  of  current  was  then  applied,  and  the  instantaneous 
resulting  Js  from  the  PDE  model  was  compared  to  the  computed 
Js  from  the  ROM.  To  simulate  the  PDE  model,  we  used  COMSOL 
Multiphysics  3.5a  [18]  coupled  with  a  MATLAB  [19]  script  to  cycle 
through  the  series  of  simulations  and  analyze  results.  Specifically, 
each  simulation  comprised  a  1  s  time  interval,  where  the  cell  cur¬ 
rent  z'app  was  modeled  as  a  Heaviside  step  function,  which  was 
applied  half-way  through  the  interval.  We  found  that  the  initial 
rest  interval  facilitated  convergence  of  the  solution  by  allowing  the 
PDE  solver  to  adjust  its  initial  conditions  before  applying  the  step 
current.  (Even  so,  convergence  of  the  PDE  simulations  proved  trou¬ 
blesome,  and  required  user  vigilance  to  ensure  reliable  results.)  The 
ROM  results  were  computed  by  a  MATLAB  script  using  Eq.  (12). 

The  cell  parameters  that  we  used  in  the  simulations  match  those 
used  in  [9,20]  and  are  listed  in  Table  1.  For  the  full-order  PDE  sim¬ 
ulations,  the  applied  current  was  varied  between  0  A  and  5.4  A  (a 
3  C  rate)  in  steps  of  0.1  A;  the  initial  cell  SOC  was  varied  between 
0%  and  100%  in  steps  of  2%,  and  temperature  was  varied  between 
-35  °C  and  45  °C  in  steps  of  20  °C.  The  solver  absolute  tolerance  was 
set  to  10-5  and  the  relative  tolerance  to  10-6.  The  default  direct 
UMFPACK  solver  was  used,  with  280  mesh  points  in  the  ID  battery 


model  and  10  624  elements  in  the  2D  electrode  model.  The  solver 
timestep  was  set  to  5  ms.  A  total  of  14025  simulations  were  run. 

For  the  reduced-order  simulations,  which  run  much  more 
quickly,  the  applied  current  was  varied  between  0A  and  5.4  A  in 
steps  of  0.05  A;  the  initial  cell  SOC  was  varied  between  0%  and  1 00% 
in  steps  of  1%,  and  the  temperature  was  varied  between  -35  °C 
and  45  °C  in  steps  of  10°C.  A  total  of  112  200  simulations  were 
run.  As  one  point  of  comparison,  the  set  of  14025  full-order  PDE 
simulations  took  more  than  eight  days  to  complete  on  an  Intel  i7 
processor,  while  the  set  of  112  200  ROM  simulations  took  a  total 
of  about  2.6  s  to  complete  on  the  same  machine.  The  speedup,  on  a 
per-simulation  basis,  is  more  than  2  000  000:1.  This  is  the  primary 
advantage  of  the  ROM  over  the  PDE  model. 

Fig.  1(a)  shows  room-temperature  side-reaction  current  density 
Js  as  computed  by  the  reduced-order  model  (which  we  now  denote 
as  Js, rom  )•  Fig-  1(b)  shows  a  compilation  of JSi  rom  over  a  range  of  tem¬ 
peratures.  We  see  two  trends  that  match  experience:  degradation 
is  worst  at  high  SOC  and  high  charge  rates. 

Fig.  2  shows  results  of  one  PDE  simulation.  This  example  was 
conducted  at  25  °C,  50%  SOC,  and  by  applying  a  1  C  charge  pulse 
at  t  =  0.5  s.  The  figure  shows  the  raw  output  of  the  simulation,  as 
compared  to  the  ROM.  Both  the  PDE  and  ROM  solutions  have  a 
non-zero  negative  side  reaction  flux  Js  even  when  the  cell  is  at  rest. 
This  is  due  to  the  fact  that  we  have  relaxed  Assumption  2  of  the  SEI 
growth  model  in  Section  2  to  also  allow  for  the  side  reaction  when 
current  in  the  external  circuit  is  zero.  Fig.  2  shows  that  the  ROM 
matches  both  the  rest  SEI  side-reaction  rate  and  the  charge-pulse 
SEI  side-reaction  rate. 

Plotted  on  the  same  scale,  the  full  PDE  solution  results  are  indis¬ 
tinguishable  from  the  ROM  results.  So,  for  comparison  purposes,  we 
define  a  relative  error  between  the  results  as 


J  is, PDE  —  is, ROM  inn 

is, err  %  =  - j -  *  100, 

is,  PDE 


where  is,pDE  is  chosen  to  be  the  value  of  Js  from  the  PDE  solution 
immediately  after  the  application  of  the  current  pulse.  Fig.  3  plots 
the  relative  error  between  the  PDE  and  ROM  solutions  for  all  25  °C 
simulations.  Between  10%  and  90%  SOC,  typical  extremum  oper¬ 
ating  conditions  for  electric  vehicle  battery  cells,  for  example,  the 
maximum  relative  error  was  0.44%. 

To  further  illustrate  the  performance  of  the  ROM  and  to  see  the 
dependence  of  SEI  layer  growth  rate  on  SOC  and  charge  rate,  Fig.  4 
plots  these  results  in  a  different  format  (25  °C).  Frame  (a)  shows 
the  side-reaction  current  density  as  a  function  of  SOC  at  different 
charge  rates  (lines  plotted  from  0C  to  3C  in  steps  of  0.5  C),  and 
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Instantaneous  degradation  rate,  25°C  Temperature  dependent  degradation  rate 


(a)  (b) 

Fig.  1.  Reduced-order  model  results  for  Js:  (a)  at  25  °C;  (b)  at  various  temperatures. 


Comparing  ROM  to  PDE  solutions 


Fig.  2.  An  example  PDE  simulation  result. 


Frame  (b)  shows  the  side-reaction  current  density  as  a  function  of 
charge  rate  at  different  SOCs  (lines  plotted  from  0%  SOC  to  100% 
SOC  in  steps  of  10%  SOC).  In  all  plots,  the  corrected  PDE  result  is 
drawn  as  a  solid  line,  and  the  ROM  result  is  drawn  as  a  dashed  line. 
In  most  cases,  it  is  impossible  to  visually  distinguish  between  the 
PDE  and  ROM  results. 

Fig.  5  shows  additional  effects  on  relative  error.  Frame  (a)  shows 
how  relative  error  varies  with  temperature.  We  see  that  the  ROM 
predictions  are  best  at  high  temperatures,  and  less  good  at  low 
temperatures.  Worst-case  Js,e rr%  in  the  10-90%  SOC  range  varies 
from  0.41%  at  45  °C  to  0.55%  at  -35  °C.  Frame  (b)  investigates  the 
effect  of  At  on  the  results.  Instead  of  selecting  the  value  for  Js,pde 
immediately  after  the  application  of  the  current  pulse  JsPDE  is  now 
selected  to  be  the  PDE  solution  0.5  s  after  the  application  of  the  cur¬ 
rent  pulse,  at  the  t  =  1  s  point.  We  see  from  Fig.  2  that  we  expect  a 
somewhat  different  result  using  this  method;  indeed,  Fig.  5(b)  ver¬ 
ifies  this  observation.  The  relative  error  is  once  again  worst  at  low 


Relative  error 


Fig.  3.  Relative  errors  between  PDE  and  ROM  solutions. 


temperatures  and  low  values  of  SOC  (where  the  absolute  amount  of 
degradation  is  small).  Relative  errors  over  1 0%  are  observed  in  some 
cases,  but  in  the  ranges  of  SOC  most  important  for  control,  where 
SOC  is  greater  than  25%,  the  worst-case  Js,err%  is  far  less,  varying 
from  0.85%  at  45  °C  to  1 .04%  at  -35  °C. 

Fig.  6  investigates  the  effect  of  a  prolonged  constant-current 
charge  at  a  1  C  rate,  as  might  be  experienced  when  a  cell  is  being 
charged.  The  PDE  is  simulated  for  3000  s,  starting  with  the  cell  at 
rest  at  1 0%  SOC,  and  1 D  profiles  of  Js{x)  across  the  anode  are  plotted 
at  time  steps  100  s,  1000  s,  2000  s,  and  3000  s.  Overlaid  on  the  plot 
are  the  average Js  values  predicted  by  the  ROM  at  that  SOC  level,  and 
the  actual  averaged Js  values  (averaged  over  the  ID  electrode)  from 
the  PDE  solution.  In  the  ROM  simulation,  the  SOC  is  updated  on  a 
second-by-second  basis  to  achieve  the  present  SOC  at  every  point, 
which  is  used  to  compute  the  value  of ]s  using  the  method  explained 


Degradation  as  a  function  of  rate 


Fig.  4.  Degradation  as  functions  of  SOC  and  rate. 
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Effect  of  temperature  on  results  Effect  of  0.5s  At  on  results 


(a)  (b) 

Fig.  5.  Additional  effects  on  relative  error:  (a)  various  temperatures:  (b)  a  different  At. 
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Appendix  A.  Model  equations 

The  following  is  a  summary  of  the  governing  equations  for  the 
model  developed  in  [4]  and  used  in  [9]  for  ideal  cell  dynamics.  The 
conservation  of  lithium  in  a  single  particle  is  modeled  by  Fields  law 
of  diffusion 


0  0.2  0.4  0.6  0.8  1 

Dimensionless  distance  across  anode 


dcs  _  Ds  d  f  2  dcs  \ 
dt  r2  dr  \  dr  J  9 


^  -°-2 
CO 
I 

E 

°  -0.4 

< 

“  -0.6 


-0.8 


100  s 


1000  s 


2000  s 


3000  s 


1 D  PDE  computation  of  Jg 

_  _  .  Avg.  over  width  of  PDE  computation  of  Jg 
,  _ ,  _  ROM  computation  of  average  J 


Fig.  6.  Comparing  PDE  and  ROM  solutions  over  long  1  C  charge. 


in  this  paper.  We  see  that  even  over  prolonged  constant-current 
charge  profiles,  the  ROM  is  accurate,  indicating  that  Assumption  1 
of  ROM  (cf.  Section  3)  is  a  reasonable  assumption  to  make. 

5.  Conclusion 

From  the  results  above,  there  is  a  clear  advantage  to  having  a 
reduced  order  model.  The  reduced  order  model  takes  the  very  dif¬ 
ficult  and  time  consuming  procedure  to  calculate  the  governing 
equations  necessary  to  obtain  Js  in  a  PDE  model,  and  reduces  them 
to  a  very  fast  algebraic  expression,  Eq.  (12),  that  a  small  micro¬ 
controller  can  easily  calculate.  As  noted,  the  speedup  is  more  than 
2  000  000 : 1  when  using  the  ROM  versus  the  PDE  solver.  Further, 
from  this  simple  calculation  of  Js,  an  incremental  model  of  degra¬ 
dation  is  obtained,  Eqs.  (13)  and  (14).  The  results  themselves  are 
very  accurate,  resulting  in  about  0.5%  relative  error  at  the  instant  of 
applied  current  to  about  1%  relative  error  0.5  s  later,  for  SOC  ranges 
of  interest. 

Due  to  the  high-fidelity  approximation  of  Js,  we  expect  that  con¬ 
trol  actions  taken  using  JsROm  will  closely  approximate  those  taken 
usingJSipDE.  Therefore,  this  model  shows  potential  for  being  used  in 
an  optimal  charger  design.  This  would  allow  a  battery  to  be  charged 
very  quickly,  but  also  allow  the  microcontroller  to  optimize  the 
charging  profile  easily  in  order  to  avoid  unnecessary  degradation. 
In  a  future  paper,  the  authors  will  explain  their  current  research  into 
optimal  profiles  that  consider  SEI  growth  as  well  as  other  degrada¬ 
tion  mechanisms. 
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with  boundary  conditions 

-Ji 


where  cs  is  the  solid-phase  lithium  concentration  (the  subscript  s 
denotes  the  solid  phase),  Ds  is  the  solid  phase  diffusion  coefficient, 
Ji  is  the  volumetric  rate  of  electrochemical  reaction  at  the  particle 
surface,  as  is  the  specific  interfacial  surface  area,  and  F  is  Faraday’s 
constant.  For  a  spherical  particle  of  active  material  that  has  a  radius 
R  in  an  electrode  with  volume  fraction  of  active  material  £s ,  the 
interfacial  surface  area  is  as  =  3 £SIR. 

Conservation  of  lithium  in  the  electrolyte  phase  gives  the  fol¬ 
lowing  equation: 


dcs 

dr 


dr 

=  0  and  Ds^ 
dr 


d(sece ) 
dt 


d_ 

dx 


+ 


1 


which  has  the  following  zero  flux  boundary  conditions  at  the  cur¬ 
rent  collectors, 


dce 

dx 


=  0 , 


x=L 


where  ce  is  the  electrolyte  phase  lithium  concentration  (the  sub¬ 
script  e  denotes  the  electrolyte  phase),  £e  is  the  electrolyte  phase 
volume  fraction  and  is  the  transference  number  of  Li+  with 
respect  to  the  solvent  velocity,  which  is  assumed  constant.  The 
effective  diffusion  coefficient  is  calculated  from  a  reference  coeffi¬ 
cient  using  the  Bruggeman  constant,  Dff  =  Des ^rug.  The  Bruggeman 
relation  takes  into  account  the  tortuosity  of  the  path  the  Li+  ions 
take  through  the  electrode. 

Charge  conservation  in  the  solid  phase  in  the  electrodes  is 
described  by  Ohm’s  law, 


-Ji  =  0, 
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with  boundary  conditions  at  the  current  collector, 


aeff  90s 

dx 


x=0 


=  aeii  90s 

dx 


x=L 


I  _  . 

A  ~  'app’ 


where  <ps  is  the  potential  of  the  solid,  aef{  =  a+ug  is  the  effective 
conductivity  of  the  solid,  A  is  the  electrode  plate  area,  and  I  is  the 
applied  current,  where  a  positive  current  discharges  the  battery. 
There  must  be  zero  electronic  current  at  the  separator  boundary, 
so 


90s 

dx 


X=Ln,s( 


90s 

dx 


=  0. 


x=Lv, 


Conservation  of  charge  in  the  electrolyte  gives  the  equation 

+h  =  0, 


where  0e  is  the  potential  of  the  electrolyte  phase,  and  /ceff  is  the 
effective  ionic  conductivity,  calculated  from  the  Bruggeman  rela¬ 
tion, 


K 


K 


eff  =  K£herug  =  (4.1253  X  1(T4  +  5.007ce  -  4.7212  x  10 
+1.5094  x  10 6c|  -  1.6018  x  \08cf)sbemg 
eff  _  2Kg77ceff  /  o  \  /  dln(/±)\ 

F  Vr+  dln(Ce);’ 


where  Rg  is  the  universal  gas  constant,  T is  the  temperature  and  f±  is 
the  activity  coefficient,  assumed  here  to  be  constant  (its  derivative 
is  zero).  This  equation  has  boundary  conditions 


90e 

dx 


x=0 


90e 

dx 


=  0. 

x=L 


The  PDEs  are  coupled  via  the  intercalation  current  density 


Ji  =  ajioj 


exp 


^aj1 


V)  -  exp  - 


acJF 

RgT 


flj 


J  =  n,p 


which  is  driven  by  the  overpotential 

rij  =<j)s-<t)e-  ujef  -  ^+Rfilm, 


where  i o  is  the  exchange  current  density, 

ioj  =  kj  (c™*  -  csj)“-  (cSJ)“CJ'(cef  *J,  j  =  n,p 


and  Ujef  is  the  equilibrium  potential  which  is  evaluated  as  a  function 
of  the  solid  phase  concentration  at  the  surface  of  the  particle.  For 
this  work, 

UTnef  =  0.7222  +  0.13876>n  +  0.029 0^/2  -  °’°172  +  0  0019 

On  0n- 

+  0.2808  exp  (0.90  -  156>n)  -  0.7984  exp  ( 0.4465+,  -  0.4108)  (15) 

r rref  (-4.656  +  88.669^  +  342.909^  +  342.9090®  -  462.4710^  +  433.4340^°) 
p  ~  (-1  +  18.93302  -  79.53204  +  37.3110®  -  73.0830®  +  95.960;°)  ' 
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